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Abstract 

We apply a molecular approach to quantify the level of hymenopteran parasitoids infestation in the lar- 
vae of the marsh fritillary (Euphydryas aurinia), a declining butterfly species, in western Bohemia, Czech 
Republic, in two subsequent years. We used the novel primer HymR157 in combination with known 
universal 28SD1F to establish a PCR detection system which amplifies hymenopteran parasitoids, but not 
the lepidopteran host. In the 14 sampled E. aurinia colonies, the infestation rates per individuum were 
33.3% and 40.2%; whereas per sampled larval colony, these were on average 38.5% (range 0-100) and 
40.1% (0-78). The per-colony infestation rates correlated with the numbers of larval webs censused per 
colony the year prior to sampling the parasitoids, pointing to a time lag in parasitoid infestation rates. The 
levels of the hymenopteran parasitoid prevalence are thus relatively high, supporting the importance of 
parasitoids for the population dynamics of the threatened host. The detection primers we developed can 


detect a range of hymenopteran parasitoids on other butterfly hosts. 
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Introduction 


Hymenopteran parasitoids are one of the most diverse groups of animals in terrestrial 
ecosystems and play a key role in the natural regulation of their host populations (La 
Salle and Gauld 1991; Forbes et al. 2018). The impact of parasitoids on their hosts can 
vary depending on their ecological specialisation, but in general they are known to cause 
significant levels of mortality in their hosts (Hawkins 1994). High levels of parasitism 
may also pose a potential threat to many threatened butterfly species, especially to vari- 
ous specialists in fragmented landscapes (cf. Anton et al. 2007). Species of the genus Eu- 
phydryas Scudder, 1872 (Nymphalidae: Melitaeini) represent a suitable system for study- 
ing host-parasitoid interactions, as egg clutches and webs with gregarious caterpillars can 
easily be detected in the field (Stamp 1981; Hula et al. 2004; Johansson et al. 2019). 
Previous studies using rearing have shown that the level of parasitism varied annually 
and depended mainly on weather conditions, which play a key role in the synchronisa- 
tion between larval development and the emergence of parasitoids (Porter 1979, 1983). 
The Marsh Fritillary, Euphydryas aurinia (Rottemburg, 1775) is an EU-protected 
butterfly, declining in many European countries (van Swaay et al. 2010), including the 
Czech Republic (Hejda et al. 2017). It belongs to a genetically polymorphic group of 
closely related taxa, the “E. aurinia complex” (cf. Korb et al. 2016), with a wide Palearc- 
tic distribution and regional habitat and host plant specificity (e.g., Munguira et al. 
1997; Singer et al. 2002; Junker et al. 2010; Korb et al. 2016). In Central and Western 
Europe, its main habitats are oligotrophic grasslands, and the most frequently used 
host plant is Succisa pratensis Moench (Dipsacaceae) (Warren 1994; Anthes et al. 2003; 
Konvicka et al. 2003; Meister et al. 2015). The butterfly is monovoltine, with flight pe- 
riod in late spring/early summer when the mated females oviposit on leaf rosettes of the 
host plant. The larvae feed gregariously in silken webs on the plants until overwintering. 
They enter hibernation with the host plants’ senescence in mid-September and resume 
feeding solitarily in April. In the Czech Republic, the distribution is restricted to the 
western part of the country (Fig. 1), where it forms three distinct metapopulation clus- 
ters inhabiting ~90 separate oligotrophic meadow patches interconnected by dispersal 
(Zimmermann et al. 2011; Junker et al. 2021; Tajek et al. 2023). This system is moni- 
tored annually by counting larval nests (cf. Ojanen et al. 2013) and displays remarkable 
within-site and inter-annual dynamics with booms and bursts (John et al. in rev.). 
Like many other insects, E. aurinia hosts numerous hymenopteran parasitoids 
(Wahlberg et al. 2001; Eliasson and Shaw 2003; Stefanescu et al. 2009). The braconids 
of the genus Cotesia (Cameron, 1891), gregarious endoparasitoids of Lepidoptera, can 
be considered the most important and numerous. Their adult females oviposit into 
haemolymph of lepidopteran caterpillars; their larvae feed internally, break through 
the cuticle in the prepupal larval instar, and form silky external cocoons, in which they 
pupate, and from which the adult wasps hatch (Kester and Barbosa 1991; Pakarinen 
2011). It was long believed that the main hymenopteran parasitoid of European Meli- 
taeini butterflies was Cotesia melitaearum (Wilkinson, 1937), remarkable for its pluriv- 
oltine development, in which successive wasp broods oviposit on successive caterpillar 
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Figure |. The distribution of Euphydryas aurinia in the Czech Republic, historic records included, based 
on Benes et al. (2002), with actualisations. The inset in upper right corner shows the position of the 
country in Europe. 


instars (Shaw et al. 2009; Pakarinen 2011). A molecular approach complicated the mat- 
ter by revealing that C. melitaearum is a complex of several cryptic species (Kankare et 
al. 2005c; Stefanescu et al. 2009). Regardless, the following hymenopteran parasitoids 
of E. aurinia have been so far reported from the Czech Republic: Braconidae — Cote- 
sia melitaearum (Wilkinson, 1937), C. tibialis (Curtis, 1830); Pteromalidae — Ptero- 
malus puparum (Linnaeus, 1758); Ichneumonidae — [chneumon emancipatus (Wesmael, 
1845), Ichneumon gracilicornis (Gravenhorst, 1849) (Shaw et al. 2009; Yu et al. 2012). 

Diverse methods were used so far to study the parasitoids of E. aurinia, and re- 
lated butterflies, ranging from field counts of hymenopteran cocoons (Ford and Ford 
1930; Porter 1983), captive rearing (Eliasson and Shaw 2003), field experiments with 
captive-reared material (Stamp 1981) to population genetic studies targeting parasi- 
toid adults (Lei and Hanski 1997; Van Nouhuys and Lei 2004). However, the ques- 
tion pivotal to the butterfly population dynamics and conservation, that of infestation 
rates relative to population cycle and state of the butterfly colonies, seems to be little 
explored. This is probably due to the work requirements for rearing both butterflies 
and parasitoids (cf. Klapwijk and Lewis 2014), combined with destructivity of such 
methods for field populations. To quantify the parasitism rates in the Czech Republic 
populations of E. aurinia, we developed a molecular method, allowing rapid and low- 
cost detection of Hymenoptera parasitoids’ incidence. 
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DNA-based methods are increasingly used in studies of parasitoid-hosts interac- 
tions (Zhu et al. 2019; Jeffs et al. 2021). The protocols so far developed for Lepi- 
doptera/Hymenoptera systems mainly focused on COI locus, commonly referred as 
barcode (Folmer et al. 1994; Hebert et al. 2004), with the hope that the solid barcod- 
ing databases will assist species* identification (Toro-Delgado et al. 2022). Particularly 
good results were obtained via a reversal approach when the adult parasitoids were 
screened for host DNA shortly after their emergence (Rougerie et al. 2011) or in a spe- 
cies-poor natural system (high Arctic: Wirta et al. 2014). However, use of COI-based 
primers may be unreliable without subsequent sequencing, because deeply phyloge- 
netically conserved bases are few and too far between in COI to place a group specific 
primer. Therefore, it was necessary to find a novel primer or primer pair which would 
amplify Hymenoptera but not Lepidoptera in the mixed samples containing the DNA 
of known lepidopteran host and unidentified hymenopteran parasitoids. We found 
such a potential primer in the nuclear region encoding the 28S ribosomal DNA. The 
novel primer, together with a primer published by Larsen (1992), targets part of the 
28S gene and aims to amplify only Hymenoptera. 

In this paper, we quantify Hymenoptera parasitoids infestation rates in a selection 
of the Czech Republic populations of Euphydryas aurinia and relate the infestation 
rates to the stage of the butterfly population cycle. Additionally, we document utility of 
our primers’ combination for rapid Hymenoptera infestation assessment in butterflies. 


Material and methods 


We sampled £. aurinia caterpillars in western Bohemia (Fig. 1) in late August and early 
September. We sampled two caterpillars per larval web (105 webs in total) from 13 
sites in 2019 and four per web (90 webs in total) from 9 sites in 2020. 

While sampling the caterpillars, we recorded the following: Julian date, to account 
for infestation changes during larval period; longest and shortest dimension of the lar- 
val web (cm); sward height, i.c., visually estimated height of surrounding vegetation in 
2.5 metre radius circles around each larval web sampled; host plant density, expressed 
as the number of Succisa flowerheads in the circle; and webs density, expressed as the 
number of larval webs in the circle. 

The material was stored in 96% ethanol, the DNA was extracted using the Tissue 
Genomic DNA Mini Kit (GenAid Biotech, Taiwan) following the manufacturer's protocol. 

We targeted a part of the 28S D1 region. We used primer 28SD1F (GGG- 
GAGGAAAAGAAACTAAC; Larsen et al. 1992) in combination with a new primer 
HymR157 (TGGCCCCATTCAAGATGG) with a resulting product of 164-167 
bp. For the primer design, we assembled a library of target sequences (Hymenop- 
tera parasitoids) and non-target sequences (Lepidoptera) from sequences available 
in GenBank (primarily PopSet 300390962, Heraty et al. 2011). We aligned the 
sequences in GENEIOUS PRIME 2020.2.4 (https://www.geneious.com) software 
and used AMPLICON software (Jarman 2003) to identify sections with concentrat- 
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ed nucleotides that consistently differ between the target and non-target groups. We 
then manually screened these promising sections and used general rules of thumb to 
design candidate primers. We aimed for at least three differences between target and 
non-target groups in the first five positions at the 3’ end of the primer for reliable 
specificity. We then tested the candidate primer in an in-silico PCR in GENEIOUS 
and optimized melting temperature in the primer pairs by extending or shortening 
the primers. 

Each in vitro PCR reaction contained 6.5 pl of Combi PPP Mastermix (Top-Bio, 
Czech Republic), 4.5 ul of HOv0,5 ul of both reverse and forward primer, and 1 ul 
of DNA template. The cycling conditions of PCR were as follows: 94 °C of initial 
denaturation (5 mins), 30 cycles at 94 °C denaturation (40 secs), 50 °C annealing (30 
secs), and 72 °C elongation (1 min), with the final elongation at 72 °C (5 mins). The 
presence/absence of PCR products was checked using agarose electrophoresis (1.5% 
gel, 150V, 30 mins). The samples with a band on the gel were assumed as positive; i.e., 
individuals infested by parasitoids (Fig. 2). 

To test the utility of the primer used, we also carried out control reactions which 
contained DNA extracted from various adult hymenopteran parasitoids and butterflies 
(Fig. 2) to assure that we amplified only the potential parasitoids and not the butterfly. 
Some of the obtained PCR products (n=4) of parasitoids from positive E. aurinia sam- 
ples were sequenced in SEQme (Czech Republic) to confirm hymenopteran origin. We 
checked the identity of the obtained sequences by using BLAST (nBLAST algorithm) 
(https://blast.ncbi.nlm.nih.gov, Altschul et al. 1990). 

Positive results were recalculated to infestation levels per larval web (1/0 factor, in- 
fested or not) and per site (% of larvae sampled). To relate the per-web infestation level 
to larval web properties, we carried out logistic regressions (binomial error distribution; 
in R 4.2.3, R Core Team 2018) with infestation (1/0) as the dependent variable; and 
the longest and shortest web ground projections, surrounding vegetation height, Suc- 
cisa density, and number of larval webs as predictors. We used the information theory 
approach, comparing the fitted regression Akaike information criteria (AIC) with AIC 
of the null model, y-+1, and considered models with AAIC > ~2.0 as fitting the data. 

To relate the per site percentual infestation to larval counts at the sites, we used 
data from annual monitoring of the sites, ongoing since 2001 (Hula et al. 2004). Over 
this time, larval counts were obtained for roughly % of site x year combinations (John 
et al. in rev., Suppl. material 1). 


Results 


Out of 210 (year 2019) and 358 (year 2020) E. aurinia caterpillars assayed for hyme- 
nopteran DNA, we obtained 70 and 144 positive results, respectively; i.e., the total 
infestation rates were 33.3% and 40.2% per individuum. On a per-site basis, this 
translates to mean+SD / median / range 38.5£29.89 / 40 / 0-100 per cent in 2019, 
and 40.1+26.51 / 50 / 0-77.5 per cent in 2020. 
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Agrypon polyxenae 
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Cotesia melitaearum 
(Wilkinson, 1937) 
; Cotesia gades 
(Nixon, 1974) 
} Cotesia gades 

{ (Nixon, 1974) 


Ichneumonoidea 


Euphydryas qurinia 
(Rottemburg, 1775) 


Pteromalus puparum 
(Linnaeus, 1758) 


Brachymeria podagrica 
(Fabricius, 1787) 


Chalcidoidea 


Cacycerus lingeus (Stoll, [1782]) 
(Lycaenidae) 

Coenonympha tullia (Miller, 1764) 
(Nymphalidae) 

Belenois aurota (Fabricius, 1793) 
(Pieridae) 

Euphydryas maturna (Linnaeus, 1758) 
(Nymphalidae) 


Figure 2. Electrophoresis gels used to assess whether the primers used can discriminate lepidopteran hosts 
and hymenopteran parasitoids a various adult hymenopteran parasitoids (PCRs are positive) and positive (P) 
and negative (N) samples of £. aurinia b four species of adult butterflies; PCRs are negative. The adult speci- 
mens of Hymenoptera and Lepidoptera were identified by M. Rindo’, M. Konvicka, and Z. Faltynek Fric. 


The sequences of the positive samples were 124—126 bp long. The most similar 
sequences in GenBank according to nblast algorithm are those of Cotesia glomerata 
(Linnaeus, 1758), with the query identity 95.2-96.03%. The next similar sequences 
did not even reach a match of 93%. 

According to the logistic regressions, none of the recorded properties of larval webs 
were related to infestation of the web (Table 1). 

At the level of individual colonies, the infestation rates were highly variable 
(Fig. 3a). The per-site infestation levels did not correlate with larval web counts from 
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Table |. Logistic regressions relating field-measured properties of larval webs to infestation of the sam- 


pled larvae. 


Model Predictor meant+SD/median/range Coefficient Residual deviance Residual DF AIC 
Null 1° 223.2 160 225.2 
Longest dimension 4,3+2.68/4/1-13 0.038 222.9 159 226.9 
Shortest dimension 7.6£12.76/12/3—40 -0.001 146.5 159 227.2 
Null 2° 181.1 194 237.7 
Julian date 240+4.3/238/235-248 0.052 178.9 193 259.0) 
Sward height 42+23.1/4/5-100 0.010 178.8 193 239.8 
Host plant density 76£49.5/60/5—200 0.004 179.6 193 237 of 
Webs density 1.7£1.98/1/0-9 -0.095 179.6 193 236.9 


The fitted single-term models are compared with the null model(s) following the information theory approach. None 
of the models was significant, as AAIC (null — fitted model) were always < 2.0. 

YWe fitted two null models, because measurements of E. aurinia larval webs were not available for 34 webs, which were 
disintegrated at the time of sampling the caterpillars. 
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Figure 3. Per-site hymenopteran parasitoids infestation rates in colonies of the butterfly Euphydryas 
aurinia in two consecutive years (2019-20), with information of caterpillar web counts in the respective 
colonies in 2018-2021 (above), and illustration of the relationship between hymenopteran parasitoids 


infestation rates and E. aurinia caterpillar web counts in the previous year (below). 
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the same year (Spearman's r= 0.11, bara = 0.46, p = 0.65) or with web counts in the 


subsequent year (r = 0.31, bocapy = 1.32, p = 0.21), but did correlate positively with the 
larval webs’ counts from the previous year (i.e., 2018 web counts for 2019 sampling, 
and 2019 web counts for 2020 sampling: 7, = 0.46, Eagar = 2-07; p = 0.054) (Fig. 3b). 

Because the absolute values of web counts highly varied among the sites (Fig. 3a), 
depending, e.g., on the site areas, we also recalculated web counts into percentage frac- 
tions of a ten-year (201 1—2021) maximum for the given site, and recalculated the cor- 
relations, with no significant results (identical year percentage web counts: r. = —0.23, 
barat = —0.99, p = 0.37; previous year percentage web counts: r= Ook, Lagig- = lads 


p = 0.21, subsequent year percentage web counts: r = 0.10, bosap = 0.40, p = 0.69). 


Discussion 


The novel combination of primers 28SD1F (Larsen 1992) and HymR157 allowed us 
labour-efficient detection of high infestation rate in the declining EF. aurinia butterfly 
by Hymenoptera parasitoids. Sequencing a selection of the obtained PCR products 
suggested that some of the parasitoids belong to the genus Cotesia Cameron, 1891. 
This is supported by the fact that other hymenopteran parasitoids known from our 
region, Ichneumon spp. and Pteromalus spp., attack pre-pupation larvae and pupae, re- 
spectively, whereas we worked with pre-diapause larvae. The gene 28S D1 is currently 
little represented in nucleotide databases for genus Cotesia, which, together with the 
unsettled taxonomy of Cotesia wasps infecting Melitaeinae butterflies (cf. Kankare et 
al. 2005a, b), precludes specific identification at this moment. This highlights the need 
to continue building comprehensive reference libraries for species identification. Such 
libraries are currently well developed for the COI gene, but supplementing it with a 
marker in another locus could improve discrimination power in some complex cases. 
Still, our approach allowed quantifying hymenopteran infestation rates in the declin- 
ing butterfly, revealing that per-colony infestation rate is affected by larval webs density 
in the previous season. 

Although the role of parasitoids on population fluctuations of FE. aurinia, and 
related species, had been proposed almost a century ago (Ford and Ford 1930), rel- 
atively few authors quantified the natural infestation rates, with widely varying re- 
sults. Infestation rates of <5% were reported, e.g., for the American congeneric species 
E. editha (Boisduval, 1852) (Singer and Erlich 1979) and E. chalcedona (Doubleday, 
1847) (Lincoln et al. 1982). Similar results were obtained for £. aurinia from Sweden, 
with rates 2.6% (Eliasson and Shaw 2003), and Spain, where the rates were 2.4-5.1% 
(Stefanescu et al. 2009). The latter study in fact covered a newly recognised species, 
E. beckeri (Lederer, 1853), feeding on Lonicera spp. The rates detected by us, 33.3% 
and 40.2% in two consecutive years, are more comparable to the situations report- 
ed for American E. phaeton (Drury, 1773) (up to ~10% in larvae prior to diapause) 
(Stamp 1981), E. maturna (Linnaeus, 1758) in the Czech Republic (69%) (Dolek et 
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al. 2006) and Sweden (32%) (Eliasson and Shaw 2003). For E. aurinia, values higher 
than ours (90%) were reported by Ford and Ford (1930) from Britain during peaks of 
cyclic fluctuation of the butterfly population, whereas Klapwijk and Lewis (2014) re- 
ported a high range of infestation rates within individual webs (4—83%) from Britain. 
This all points to a high variation among sites, seasons, parasitoids’ generations, and 
Euphydryas species in hymenopteran infestation rates. In detailed studies of the closely 
related Melitaeinae model species, Melitaea cinxia (Linnaeus, 1758), this variation was 
attributed to spatial positions of butterfly colonies, competition among parasitoids 
and hyperparasitism, and annual variation in weather (e.g., Lei and Hanski 1997; Van 
Nouhuys and Lei 2004). Arguably, some of the reported variation may also be due 
to the diversity of methods applied by various authors, ranging from field counts of 
infested caterpillars (Lei and Hanski 1997), through captive rearing (Stamp 1981; Eli- 
asson and Shaw 2003; Klapwijk and Lewis 2014), to molecular methods as used here. 
Possibly, the molecular detection reveals higher infestation rates than rearing, because 
some of the infested larvae may die prior to the parasitoid emergence. Causes of this 
mortality are then interpreted as “unknown” (e.g., in Eliasson and Shaw 2003). 

Klapwijk and Lewis (2014) observed that the probability of caterpillar web infesta- 
tion increased in webs isolated from other F. aurinia larval webs. We did not detect any 
relationship to the webs or surrounding vegetation parameters, but these results may be 
biased, because — for conservation concerns — we sampled the caterpillars solely from 
colonies containing a high number of larval webs during the sampling. The mean+SD 
webs’ number of sampled colonies were 116£97.3 and 6066.6, whereas the numbers 
across all colonies were 53+57.8 (n = 44) and 19+38.0 (n = 70) in 2019 and 2020, 
respectively (John et al. in rev.). Conservation concerns also prevented quantifying the 
proportion of infestation per larval web, which would require killing all the caterpillars 
(cf. Klapwijk and Lewis 2014). 

It also should be noted that our approach did not distinguish parasitoids from 
hyperparasitoids; i.e., the insects developing within parasitoids and thus killing them 
(Nair et al. 2016). A hyperparasitoid, however, can infest only a larva already infested 
by a parasitoid, and hence hyperparasitoids presence does not affect our findings on 
hymenopteran infestation rates. 

With all the limitations, we found that the infestation rate per site positively corre- 
lated with per-site caterpillar webs’ numbers of the previous year. This is fully expected 
if the parasitoids need a rich resource supply (i.e., high host density) to multiply in 
a butterfly colony, depleting the hosts’ numbers in the process (Ford and Ford 1930; 
Frazer 1954; Porter 1981, 1983). A time delay in parasitoids infestation, and higher 
likelihood of infestation of larger and more connected host colonies, were found by 
Lei and Hanski (1997) in the metapopulation system of M. cinxia and its parasitoids. 
Although the inter-annual abundance changes of Melitaeini butterflies’ colonies are 
likely influenced by numerous other factors, including variation in weather (Brunbjerg 
et al. 2017) or site vegetation management (Johansson et al. 2019; Tajek et al. 2023), 
natural enemies’ pressure certainly plays a significant role. 
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Conclusions 


The DNA-based method detected high hymenopteran parasitoids’ infestation rates in 
colonies of the declining butterfly, Euphydryas aurinia. These rates, however, widely 
varied among the butterfly colonies and between two study years, likely interfering 
with, and possibly driving, the inter-annual butterfly abundance changes within colo- 
nies, as well as the metapopulation dynamics of the butterfly, described in detail by 
John et al. (in review). Our primer combination seems to be promising for wider 
use in detecting infestation of butterflies by Hymenoptera parasitoids. In our control 
tests, it amplified various Hymenoptera but not several butterfly species tested (Fig. 2). 
However, potential users should first test that the detection works also in their system 
to avoid false positive and false negative results. 
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Explanation note: Spreadsheet with information on all 97 E. aurinia colonies moni- 
tored in western Czech Republic in 2001-2021, with caterpillar webs counts for 
individual years. The fourteen colonies sampled for hymenopteran parasitoids are 
indicated in red. 

Copyright notice: This dataset is made available under the Open Database License 
(http://opendatacommons.org/licenses/odbl/1.0/). The Open Database License 
(ODDbL) is a license agreement intended to allow users to freely share, modify, and 
use this Dataset while maintaining this same freedom for others, provided that the 
original source and author(s) are credited. 
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